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ANALYSIS OF THE TRANSIENT FLOW OF HYDROGEN IN HEATED 
PASSAGES IN THE LAMINAR -TURBULENT TRANSITION REGION 

by Robert W. Leko 
Lewis Research Center 

SUMMARY 

A transient analysis has been conducted of the flow characteristics in a single heated 
tube of circular cross section. The study was made in association with phenomena occur- 
ring in the nuclear rocket reactor core. Twenty-one cases were considered involving the 
following: twelve perturbations of heat flux into the tube wall, five perturbations in inlet 
pressure, and four perturbations of pressure drop. The results of the analysis are de- 
picted in terms of the phase plane that generalizes the results and illustrates all flow 
behavior of the laminar instability condition. 

The following four external parameters, which influence the transient behavior of the 
system, were considered: inlet pressure, outlet pressure, inlet enthalpy, and heat flux 
into the tube wall. 

The analysis has shown that there are five kinds of transients. The first type of 
transient corresponds to an initial laminar flow rate on the steady -state curve which 
when perturbed goes to a turbulent equilibrium point. The second type of transient starts 
from the laminar flow rate on the steady -state curves and terminates at zero flow rate. 

A third sort of transient begins from an operating point corresponding to a turbulent flow 
rate on the steady -state curve and proceeds to an equilibrium point at a higher turbulent 
flow rate. The fourth kind of transient is initiated from a turbulent point on the steady- 
state curves and reaches a stable equilibrium point at a lower turbulent flow rate. 

Finally, the last kind of transient begins at a turbulent steady -state point and terminates 
at zero flow. These various transients, which are approximately independent of past 
initial conditions, can be more simply described by a family of curves in a phase plane 
plot. 

The "stability criterion" near equilibrium points is simply that the four external 
parameters give an operating point with a positive pressure drop as a function of flow 
rate slope on the characteristic curves. If this slope is negative or if there is no operat- 
ing point, then the condition is unstable. 



INTRODUCTION 


The shutdown transient for the nuclear rocket after a period of full power operation 
lasts for a period of several days to remove the decay heat. Hydrogen propellant is used 
as a coolant. It is essential to use a minimum propellant for this purpose in order to in- 
crease the useful payload. When the minimum flow rate of coolant consistent with the 
maximum safe material temperature limits of the reactor core is used, a possible flow 
problem may result. 

The core of a nuclear rocket consists of thousands of closely spaced parallel pas- 
sages. To simulate the problem that may arise in the nuclear rocket, it has become 
conventional practice to use a single tube with a plenum chamber at each end as a model. 
It has been demonstrated analytically in reference 1 that in the laminar -turbulent transi- 
tion region two possible mass flow rates may exist in a single heated tube for a given 
pressure drop and heat flux in the steady state. Also, the negative slope of the charac- 
teristic curve is associated with a potential instability. These two aspects define the 
problem that is referred to as the "laminar instability problem." The characteristic 
curves have been verified experimentally in tests on a single heated tube (ref. 2). 

Since it is possible to have two distinct flow rates for a given pressure drop across 
the reactor core, it is conceivable that one passage may have turbulent flow entering and 
changing to laminar flow while an adjacent passage may have laminar flow entering and 
leaving. Even if the flow is stable, a problem exists. The heat-transfer coefficients are 
different for turbulent and laminar flows so that the tube walls will be at different temper- 
atures. This condition will give rise to thermal stresses and may result in a failure of 
the reactor core. 

A steady -state analysis was conducted to define the influence of the external param- 
eters. A transient analysis was conducted to find the kinds of transients occurring and to 
reveal the general nature of the flow stability problem; the analysis was accomplished by 
perturbing heat flux into the tube wall, inlet pressure, and pressure drop. This report 
presents the results of these analytical analyses. 


ANALYSIS 

Steady-State Analysis 

The steady -state analysis was performed in order to determine the influence of the 
several parameters on the characteristic curves. The geometry of the tube consisted of 
a length of 52 inches, an inner diameter of 0. 1160 inch, and an outer diameter of 0. 1476 
inch. This analysis used a test section consisting of an 80 -percent nickel - 20 -percent 
chromium alloy, which has a density of 0. 304 pound per cubic inch. 
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Q(x) - constant 



Figure 1. - Schematic of simulated physical system. 


Figure 1 shows the physical system being considered. It consists of a single heated 
tube with a plenum chamber at each end. For the purposes of analysis, the tube length 
was considered divided into 26 stations each 2 inches long. 

The basic assumptions of the calculations include the following. First, it was as- 
sumed that axial heat transfer along the tube wall by conduction is negligible. Therefore, 
the agreement between the analytical results of this study and the real physical situation 
improves as the tube wall thickness is decreased. Second, it is assumed that the overall 
pressure drop is comprised of two elements. One element arises due to the friction drop 
while the other occurs by virtue of the momentum pressure drop: 


AP q - APp + 


AP 


M 


( 1 ) 


(Symbols are defined in appendix A. ) An iterative technique is used to calculate the over- 
all pressure drop. 

It is further assumed that all the heat entering the tube wall goes into an increase of 
the enthalpy of the fluid: 


AH = ALQ/W 


( 2 ) 


The wall temperature is computed from Newton’s law of convective heat transfer: 

T w (x) = ALQ/(hA) + T b (x) (3) 

A more detailed description of the steady-state analysis is included in appendix B. 


Transient Analysis 

The transient analysis involves a quasi-steady state method of approach. It is 
assumed that steady -state flow conditions exist at any instant of time corresponding to the 
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wall temperature distribution, inlet pressure and enthalpy, and outlet pressure. The re- 
sulting flow and fluid state throughout the tube results in a heat-transfer rate into the fluid 
which gives a new wall temperature distribution at the next instant of time. The basic 
assumptions of the transient analysis are the same as those of the steady-state analysis. 

Pressure drop is obtained by an iterative technique from the basic equation 

iP 0 = + aP M « 

After the pressure drop is computed for a station, the heat-transfer problem is solved 
for that station. It is then assumed that the difference between the heat flux entering the 
tube wall and that transferred to the fluid goes into an increase of temperature of the tube 
wall: 


Q^x) = Q(x) - Qfjte) (4) 

Finally, the new wall temperature distribution is computed: 

V new« ' T w, oldW + ALQwW8/(p m V m c p; m N) (5) 

This procedure is continued until the entire tube length is completed. 

A more detailed explanation of the transient analysis is given in appendix C . A digi- 
tal computer program served as a tool to perform the actual numerical computations. 


Correlations 


The friction factor and heat -transfer correlations common to both analyses are the 
following: 

Friction factors: 


f £ = 16/Re b 


f t = 0. 0014 + 0. 125 mej* 32 


Heat-transfer coefficients: 
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0.525 


h. = 0.021 - Re?* 8 
t D b 




RESULTS AND DISCUSSION 

The steady-state parametric study and the transient analysis were both performed 
with a composition of 100 percent parahydrogen as the working fluid. 

The results of the steady -state analysis are shown in figures 2 to 4. The influence 
on the characteristic curves of varying the heat flux into the tube wall while holding inlet 
pressure and enthalpy constant at the respective values of 10. 0 pounds per square inch 
absolute and 351. 0 Btu per pound is shown in figure 2. The influence on the characteris- 
tic curves of varying the inlet enthalpy while holding the wall heat flux constant at 
0. 0016 Btu per second per inch and the inlet pressure constant at 10 pounds per square 
inch absolute is shown in figure 3. The influence on the characteristic curves of varying 
the inlet pressure while holding wall heat flux and inlet enthalpy constant at 0. 0016 Btu per 
second per inch and 351. 0 Btu per pound is shown in figure 4. The single curve of fig- 
ure 4 shows that both inlet and outlet pressure effects are included in the parameter 
Pi n AP. The Pj n AP term was chosen as an ordinate for figures 2, 3, and 4 because it 
brought together the characteristic curves for various inlet pressures and thereby sim- 
plified that part of the transient analysis. 

The characteristic curves are shown to demonstrate the same shape for all param- 
eters: wall heat flux, inlet pressure, and inlet enthalpy. Because of this, it is possible 
to obtain two distinct flow rates for a given pressure drop. 

The characteristic curves are shown for flow rates from 2xl0 -8 to 2x10”^ pound per 
second. The range of interest in the nuclear rocket application is easily contained within 
this span. The lower limit of the characteristic curves is determined by wall material 
temperature limitations of the reactor core. The wall heat flux varied from 0. 0012 to 
0. 002 Btu per second per inch, which is a typical range of values encountered during the 
shutdown of the nuclear rocket. The upper limit of flow rate of the characteristic curves 
for values greater than 2X10 pound per second was determined by the sonic velocity 
limit of the working fluid. 

The three classes of problems subjected to the transient analysis were defined as 
follows: 

Class 1 - The effect of perturbations in wall heat flux 

Class 2 - The effect of perturbations in inlet pressure 

Class 3 - The influence of perturbations in pressure drop 
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(a) Increase in perturbation. 



Time, sec 

(b) Decrease in perturbation. 

Figure 5. - Form of perturbations in inlet pressure, wall heat flux, 
and required pressure drop. 


Figure 5 shows the two types of perturbations used in the analysis. Perturbations in 
wall heat flux and inlet pressure consisted of a pulse of infinite width. Perturbations of 
required pressure drop consisted of a pulse of finite width. The external wall heat flux 
distribution for both the steady-state and transient analysis was constant. 

The results of the transient analysis for all cases are presented in terms of the 
phase plane. In general, such a plot consists of a time derivative of a variable against 
the variable (ref. 3). The phase plane representation was selected since it provided a 
convenient method of presenting and summarizing on a single plot the results of several 
transient cases. The variable common to all stations, flow rate and its time derivative, 
is used. 

A typical flow rate as a function of time response is shown plotted in figure 6 for the 
first class of perturbations. The accuracy of the curve is determined by the size of the 
multiplier of the time constant specified. A multiplier of one was used throughout this 
analysis. 

The phase plane figures were generated by graphically differentiating curves of the 
type shown in figure 6. The accuracy of this method was checked by analytically deter- 
mining the derivatives. This was accomplished in the following manner. A least-squares 
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curve fit was used to approximate the data of a flow rate as a function of time response. 

A hyperbolic tangent whose argument was time was the principal function. From the hy- 
perbolic term a corrective term consisting of a cubic polynominal in time was subtracted. 
This function fit the data well except at time zero where a maximum error of 6 percent 
occurred. The result was differentiated with respect to time and the graphical analytical 
development of the phase planes compared. Satisfactory enough agreement was shown 
between the two methods to justify the graphical technique. 

Twelve cases of perturbations involving wall heat flux were studied (class 1), and 
the results are shown in figures 7 and 8. The phase plane of figures 7(a) and (b) was de- 
termined for the restriction of a constant pressure drop of AP = 0. 13 pound per square 
inch absolute. The pressure drop was held constant to simulate the process occurring in 
the reactor core of the nuclear rocket. The curves are parametric with the final value 
of external wall heat flux. 

The interrelatedness of the steady -state and transient analyses becomes apparent at 
this point of the discussion. The effect of perturbing a steady-state operating point in 
figure 2 (p. 6) by a perturbation involving a decrease in wall heat flux while holding the 
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(a) Transients originate from laminar initial mass flow rate. Initial conditions: heat flux into tube wall, 0.0016 Btu per 
second per inch; 1. 345x10 ~ 5 pound mass per second. 
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(b) Transients originate from turbulent initial mass flow rate. Initial conditions: heat flux into tube wall, 0. 0016 Btu per 
second per inch; flow rate, 3. 7849x10"^ pound mass per second. 

Figure 7. - Phase plane representation of transient analysis for wall heat flux perturbations. Constant inlet pressure and 
enthalpy of 10 psia and 351.0 Btu per pound, respectively^ pressure drop, 0. 130 psia. 



Flow rate, W, lb mass/sec 


(c) Transients originate from various mass flow rates. Initial conditions: heat flux into tube wall, 0. 0012 Btu per second 
per inch. 


Figure 7. - Concluded. 


pressure drop constant is considered. A new stable point will be reached on the turbu- 
lent branch where the final wall heat flux curve intersects the constant pressure drop 
line. This transient is shown as curve A in figure 7(a). 

An increasing perturbation in wall heat flux from the initial condition on the laminar 
branch was unstable in all cases. This resulted because the line of constant pressure 
drop does not intersect the perturbed wall heat flux curve in figure 2. The flow rate went 
to zero as shown by transient B in figure 7(a). 

A perturbation representing a decrease in wall heat flux from the turbulent initial 
condition in figure 2 reached a stable solution at a higher turbulent flow. The transient 
is shown as curve C on figure 7(b). 

One of two possibilities exists if the perturbation represents an increase of wall heat 
flux for the turbulent initial condition. If the perturbation is small enough so that the 
line of constant pressure drop is intersected by the perturbed wall heat flux curve in 
figure 2, a stable operating point can be reached on the turbulent branch at a lower flow 
rate. This type of transient is shown as curve D on figure 7(b). If the perturbation is so 
large that this condition is not met, the system drives to zero flow rate, that is, is 
unstable. This condition is shown by transient E in figure 7(b). Indeed, there are five 
kinds of transients as shown by figures 7(a) and (b) collectively. Two of these transients 
correspond to laminar initial conditions, while the remaining three correspond to turbu- 
lent initial conditions of flow rate. 
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Figure 7(c) shows that, corresponding to the four external parameters, a unique 
curve exists in the phase plane. A specific case in point is the curves in figure 7(c) 
determined by the following four external parameters: inlet pressure of 10 pounds per 
square inch absolute, outlet pressure of 9. 87 pounds per square inch absolute, inlet en- 
thalpy of 351. 0 Btu per pound, and final external wall heat flux of 0. 0012 Btu per second- 
inch. Curves for four initial conditions corresponding to flow rates of 
Wj = 0. 55990X10 -5 , W 2 = 0. 80900X10" 5 , Wg = 1. 345X10" 5 , and W 4 = 3. 7849X10" 5 pound 
per second are seen to coalesce into a single curve after the effects of initial conditions 
have fairly quickly settled out. 


11 



Figure 9. - Phase plane plot for inlet pressure perturbations. Pressure drop, 0. 13 psia ; heat flux into tube wall, 0.0016 Btu 
per second per inch; enthalpy, 351 Btu per pound mass. 


Thus, if the initial part of the transients in figure 7 are ignored, a more general 
plot, independent of initial or past conditions, can be obtained. Figure 8 is thus obtained. 
The part of the curves of figure 7 that are ignored are filled in by using the equilibrium 
conditions; that is, the zero crossings of the curves are known. This procedure is nec- 
essary because the generation of that part of the curve would be very time consuming 
using the computer code. 

If the external heat flux, inlet pressure and enthalpy, and outlet pressure are given, 
the transient will follow the corresponding curve of figure 8. An initial arbitrary wall 
temperature distribution will cause an initial deviation from those curves. This devia- 
tion is rather quickly settled out because of the strong coupling between stations along 
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Figure 10. - Phase plane plot for required pressure drop perturbations. Heat flux into tube wall, 0. 0016 Btu per second per 
inch; inlet enthalpy, 351 Btu per pound mass. 


the tube. 

The time sense on the phase plane plot is obtained by noting that the time interval to 
go from one point to another is the area under the reciprocal of the curve. Thus, because 
of the nonlinearity of the curves of figure 7 or 8, the time to go from the left equilibrium 
point to the right equilibrium point, for example, is much longer at relatively higher heat 
fluxes. 

Figure 9 is the phase plane plot for the transients generated by perturbing the inlet 
pressure. Indeed, the five cases shown represent the five kinds of transients noted 
before. The similarity of figure 9 to figures 7 and 8 is apparent. These curves can be 
considered as retracing the same kind of phase plane curves as figure 8 after initial 

condition effects have settled out. 

Figure 10 shows plots for the transients generated by pulse -type perturbations on 
outlet pressure. Each branch follows a piece of the curves of figure 8. Figure 10 shows 
how a more general perturbation (a pulse) may be traced out using figure 8. It also fur- 
ther supports the contention that the plot of figure 8 is generally true and can be used to 
predict transient behavior and stability characteristics. 
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SUMMARY OF RESULTS 


The results of the analyses described herein yielded the following principal conclu- 
sions: 

1. The nonequilibrium, nonlinear ” laminar instability” behavior can be illustrated on 
a phase plane diagram like figure 8. 

2. The curves in the phase plane showing the transient behavior are of two families: 
one family intersects the flow rate axis at two points and the origin, and the other lies 
below the axis and intersects it at only the origin. 

3. Actual transients, which can start from arbitrary wall temperature distributions, 
will approach unique curves in the phase plane corresponding to the four external param- 
eters of inlet pressure, outlet pressure, inlet enthalpy, and external wall heat flux after 
the effects of initial wall temperature distribution have settled out. 

4. Each point of the curve in the phase plane represents, approximately, a corre- 
sponding quasi-equilibrium wall temperature distribution. 

5. Five kinds of transients were obtained: 

(a) Those beginning at laminar steady-state initial condition go to turbulent equi- 
librium point. 

(b) Those beginning at laminar steady -state initial condition go to zero flow rate. 

(c) Those beginning at turbulent steady -state initial condition go to a higher tur- 
bulent flow rate equilibrium point. 

(d) Those beginning at turbulent steady -state initial condition go to lower turbu- 
lent flow rate equilibrium point. 

(e) Those beginning at turbulent steady-state initial condition go to zero flow 

rate. 

6. There were no instabilities of constantly increasing flow rate. 

7. All phase plane curves are probably similar, qualitatively, to those shown in this 
report. 

8. The stability criteria near equilibrium points is simply that the four external 
parameters give an operating point on the characteristic curve with a positive slope. If 
this slope is negative or if there is no operating point, then the condition is unstable. 

Lewis Research Center, 

National Aeronautics and Space Administration, 

Cleveland, Ohio, June 10, 1966, 

121-30-02-03-22. 


14 



APPENDIX A 


SYMBOLS 


A area, in. ^ 

Cp, m specific heat of nickel chromium 
alloy, Btu/(lb mass)(°R) 

D diameter, in. 

f friction factor 

G mass flow per area in square 

feet, lb mass/ (sec) (ft^) 

g c acceleration due to gravity, 

32. 2 (lb mass -ft)/ 

(lb force) (sec^) 

H enthalpy, Btu/lb mass 

AH change in enthalpy 

h heat-transfer coefficient, 

Btu/(sec)(in. ^)(°R) 

J mechanical equivalent of heat 

constant, 778 (ft-lb force)/Btu 

K thermal conductivity based on 

film temperature 

L length, in. 

AL station length, in. 

N multiplier of time constant 

Nu Nusselt number 

P pressure, psia 

Pr Prandtl number 

AP (inlet pressure - outlet pres- 

sure), psia pressure drop 

Q heat flux into tube wall, 

Btu/(sec)(in. ) 


Qq heat transfer to fluid, 
Btu/(sec)(in. ) 

Qw heat stored in wall, 
Btu/(sec)(in. ) 


Re 

Reynolds number 

T 

temperature 

, °R 

t 

time, sec 


m 

volume of tube wall material, in. 

W 

flow rate, lb mass/sec 

e 

time constant, 9 = minimum 


Ip V C 1 



p m m p, m 



hA 


At 

time increment A0/N, sec 

Pfl 

fluid density 

, lb mass/ft^ 

^m 

density of nickel chromium alloy, 


lb mass/in. ® 

Subscripts: 


avg 

average 


b 

bulk 


calc 

calculated 


F 

friction 


f 

final 


fi 

film 

' 

fl 

fluid 


i 

initial 


in 

inlet 


j 

3 th iteration 



15 



SL 

laminar 

o 

outlet 

M 

momentum 

P 

plenum 

m 

material (nickel chromium alloy) 

t 

turbulent 

O 

overall 

w 

wall 
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APPENDIX B 


DETAILS OF STEADY-STATE ANALYSIS 


The purpose of this analysis was to provide a computational procedure for pressure 
drop and wall temperature distribution for a single heated tube when the flow rate, geom- 
etry, inlet pressure, inlet enthalpy, and heat flux into the tube wall are known. 

The tube length is divided into discrete finite intervals called stations to facilitate 
computing the conditions that exist at a definite axial position. Since the inlet pressure 
and temperature are known, the remaining fluid properties at the inlet may be determined 
immediately. The computed Reynolds number is based on properties at the inlet to the 
station. The pressure at the outlet of a station is found by an iterative technique. For 
the first trial, the outlet pressure of the station is assumed equal to the inlet pressure of 
the station. Outlet enthalpy is calculated as 


H = H. + 
o m ^ 


(Bl) 


These values of enthalpy and pressure are used to define the state of the fluid at the 
outlet and so they define the remaining properties there. 

A new exit pressure is next calculated as 


P = P. 

o, new in 


G 2 

144 g c 



1 \ _ 4fL G 2 
p in / 144D 2gp avg 


(B2) 


If this new value of outlet pressure and the previously calculated outlet enthalpy are used, 
a new outlet state of the fluid is defined from which the remaining thermodynamic proper- 
ties can be found. 

Iteration is continued using equation (B2) until the outlet pressure converges to some 
prescribed tolerance. Once the outlet pressure has converged, it and the other proper- 
ties at the outlet of the current station become inlet conditions to the next station. Again 
the iterative technique described previously is used to find the outlet conditions for that 
station. This procedure is repeated until all the stations comprising the tube length 
have been treated. 

The fluid bulk temperature is adjusted to relate to the center of each station. Finally, 
the wall temperature distribution is computed from 

T w (x) = + T b (x) (B3) 
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where the appropriate laminar or turbulent heat-transfer coefficient is used, depending 
on the Reynolds number for that station. 

A computer program exists to perform the computational procedure described pre- 
viously. The program employed a subroutine, described in reference 4, to obtain the 
fluid properties. 
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APPENDIX C 


TRANSIENT ANALYSIS 


The purpose of this analysis was to provide a computational procedure that would 
permit flow rate, wall temperature, and other variables to be computed as a function of 
time. It was assumed that a quasi-steady-state technique could be employed. In addition, 
it was assumed that inlet pressure, inlet temperature, wall heat flux, overall pressure 
drop, geometry, and initial wall temperature distribution are known. 

The analysis treats the flow at a given instant of time assuming steady-state condi- 
tions prevail along the tube in the manner described in the following steps: 

(1) Static pressure is computed as follows (see 
sketch (a)): 



(a) 


l 


£ 


Position A Position B 

y— - j 

j-*— Station length— *-{ 

^Tl’V 


I 


P 2 , t 2 , h 2 


(b) 


2 


G 

P 1 = P -1.2 (Cl) 

V P 

(2) The static enthalpy at the entrance of the tube is 

obtained as follows: 

r 2 

H i - H p - < C2 > 

(3) The model is shown in sketch (b) for each station. 


If the station is the first station of the tube, Pj equals P calculated in step 1, Hj 
equals H calculated in step 2, and Tj is the fluid temperature in the plenum. If the 
station is not the first, properties are equal to the properties at position B from the pre- 
vious station. 

(4) The Nusselt modulus is calculated based on the conditions at position A. 

(5) The heat-transfer coefficient and heat-transfer are calculated: 

h = (C3) 

D 


Of! “ ^w ‘ T b^ 


(C4) 


(6) The enthalpy at position B is calculated: 
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I 


mu 


i in i i 


■III ■! 


AL 


H 2 - H l + 


Qfl 


W 


(C5) 


(7) It is assumed that P 2 = Pp and H 2 from ste P 6 is used to define the state at 
position B. The remaining properties at B are determined by the use of the subroutine 
program. 

(8) The pressure drop is calculated from 


P 


2 



j-Lt _2i_ 

144 ^£ c P av g 


®c \P 2 Pi J 

x y 


(C6) 


(9) Using P 9 from step 8 and H 2 from step 6 define a new state at position B 
from which T 2 and remaining properties can be obtained. 

(10) The average and T^ are calculated: 


T b, avg 


(T 2 ^ T t ) 
2 


T f i, avg 


( T w tT b,avg> 


*N 


(C7) 


(11) Go back to step 4 and base the calculations of the Reynolds and Prandtl numbers 
on the new temperatures computed in step 10. Steps 5 to 10 are repeated until two Qq's 
calculated in succession in step 5 agree within some prescribed tolerance. When this 
tolerance is reached, the station is completed. The exit pressure, temperature, and 
enthalpy become entrance properties for the next station. 

(12) The AP ca i c is calculated: 


AP calc P in 


- P, 


(13) The flow rates are then adjusted by the equation 


(C8) 


/AP 


X 


W new W old 


required \ 
AP 


n 


calc 




(C9) 


where 


20 




(14) The calculation procedures of steps 1 to 14 are repeated until ^ re q U i re( j 
agrees with AP ca j c to some prescribed tolerance. 

(15) After the pressure drops, bulk temperatures, and heat-transfer coefficients are 
obtained, the time constant for each station and new wall temperature distribution are 
obtained as follows: 


p V C 

0 _ m p, m 
hA 


T 

w, new 


Q^x) = Q(x) - Q fl (x) 


AL 0Q(x) 

M = T w old (x) + — 

w > old p V C N 
F m m p,m 


(CIO) 


(CU) 

(C12) 


A new AP based on the T w new obtained in equation (C12) is calculated by returning 
to step 1 and using the final flow rate calculated in step 13. 

(16) The calculations are continued until all T w 's change less than a prescribed 
tolerance. 

A computer program employing the mechanics described serves as a tool to perform 
the numerical calculations. 
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